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ABSTRACT* Digital filters can be realized with remarkably email computational 
burden as polyphase recursive all-pass networks. These networks are formed as a 
sum of M parallel all-pass subf liters with phase shifts selected to add construc- 
tively in the passband and destructively in the stopband. The design technique we 
present here starts with prototype M-path recursive all-pass filters obtained by 
a new optimizing algorithm (described in detail elsewhere) - These filters are 
operated with combinations of all-pass transformations, resampling, and cascading 
options. We demonstrate a number of designs using the new algorithm along with 
examples of filtering systems obtained by the techniques presented. 



1. INTRODUCTION 

A standard finite impulse response (FIR) 
filter can be modeled, as indicated in Fig. 
1, as the weighted summation of the contents 
of a tapped delay line. This, summation is 
shown in (1) and the transfer function of 
this filter is shown in (2) • 
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FIGURE 1, TAPPED DELAY LINE FIR FILTER 

We note the frequency selective characteris- 
tics of the filter output is due to the 
phase shift between successive samples in 
the tapped delay line. This phaser summation 
is shown in Fig. 2 for two different fre- 
quencies. 

Note that it is the phase shift, and not the 
weighting terms which is responsible for the 
frequency dependent gain of the filter. The 
weighting terms are used primarily to con- 
trol passband width and stopband attenuation 
of the filter. An idea derived from spatial 
beamforming (known as beam spoiling) [1) and 
from signal synthesis (peak-to-rms control) 
[2) is to use structured* phase shifts (with 
arbitrary amplitudes) to accomplish the same 
goal . 



FIGURE 2. PHASOR SUMMATION OF FIR FILTER 
OUTPUT FOR TWO FREQUENCIES 

The class of filters we describe here re- 
places the weights of the FIR filter with 
all-pass subfilters which exhibit unity gain 
with frequency dependent phase shift. In 
anticipation of the polyphase structure, the 
all-pass subfilters are, as shown in (3) and 

^ Fig. 3, first order polynomials in Z" M where 

. « is the number of taps in the structure. 
This composite form is indicated in Fig* 4a 
with a modified form reflecting a polyphase 
structure in 4b. The transfer function of 

: this structure is shown in (4). 
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FIGURE 3, H-TH ORDER ALL- PASS BUBFILTER 
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Fig. 7" show examples of equal ripple designs 
obtained from this new algorithm for a 2- 
path and for a 5-path filter with two and 
three all-pass stages per path respectively. 
Interacting constraints on the optimal de- 
sign (described in the aforementioned paper) 
limit the possible number of stages per path 
for the 5-path filter to the sequences 
(1,1,1,1,0), (2,2,2,1 { ,1) , (3,3,2,2,2), etc. , 
so that the optimum 5-path filter uses three 
less stages than allocated to the design. 
Similarly the 2-path filter is limited to 
the sequences (1,0), (1,1), (2,1), (2,2), 
(2,3), (3,3), etc. 

2. TWO PATH FILTERB 

We find a wealth of interesting properties 
and clever applications for this structure 
even if we limit our discussion to a 2-path 
filter. As designed, the 2-path filter is a 
half bandwidth filter with the 3-dB band- 
width at 0.25 f f . If the zeros are restricted 
to the half sampling frequency, the filter 
identical to a half bandwidth Butterworth 
obtained by the standard warped bilinear 
transform. For this special case, the real 
parts of the root locations go to zero. If 
the zeros are optimized for equal ripple 
stop band behavior the filter becomes a 
constrained Elliptic filter. The constraints 
are related to a property of complementary 
all-pass filters. We define the all-pass 
sections for the 2-path filter as H 0 (9) and 
^(8) and, as indicated in Fig. 8, the scaled 
by one-half sum and difference of these 
paths by A(6) and B(B) , the lowpass and 
highpass filters, respectively. We know that 
the all-pass sections satisfy (5) from which 
we derive the power relationship between the 
lowpass and highpass filters shown in (6). 



For small e., we can ignore (eA 2 which leads 
to (8). 
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FIGURE 8. COMPLEMENTARY LOWPASS AND HIGHPASS 
FILTERS FROM 2-PATH ALL- PASS NETWORK 

iH 0 (9)| 2 » 1*1,(6) [ 2 * 1 (5) 

|A(e)|*+|B(e)|* - |o.5(H 0 (e)+H t <em* (6 a) 

+ jO.5{H 0 (9)-R 1 (e))| 2 

-O.5[|H O {0) | 2 +|H,(B) J 2 ] (6b) 

- 1 (6c) 

Now to interpret how this relationship 
impacts the 2-path filter. For the comple- 
mentary filters we define the minimum pass- 
band gain of by 1-c. and the peak stopband 
gain by € ? . Substituting these gains in (6) 
we obtain (7) . 



(l-c t ) 2 + (c 2 ) 2 = 1 
(l-€ t ) ? - 1 - (f})* 
1 - 2e x + (e,) 2 = 1 - {e z ) 2 



(7a) 
(7b) 
(7c) 



f, - 0.5 (€ Z ) 2 



(8) 



Thus if the stopband attenuation is selected 
to be 0.001 (60.0 dB) , the passband ripple 
is 0.000 000 5 (126.0 dB) or the passband 
minimum gain is 0.999 999 5 (-0.000 000 22 
dB) . As can be seen these filters have a 
remarkably flat passband. We conclude that 
the 2-path equal ripple stopband filters are 
elliptic filters with coupled passband and 
stopband ripple with -3-0 dB gain at 0.25 f $ . 
There are two significant differences be- 
tween the 2-path and the standard implemen- 
tation of the equivalent elliptic filter. 
The first is a four-to-one savings in multi- 
plications; a fifth order half bandwidth 2- 
path filter requires only two coefficients 
as opposed to the direct implementation 
which requires eight (scaling included). The 
second is that all-pass structures exhibit 
unity gain to internal states hence do not 
require extended precision registers to 
store internal states as do standard imple- 
mentations. 

As with any filter design, the transition 
bandwidth can be traded for out of band 
attenuation for a fixed order or transition 
bandwidth can be steepened for a fixed 
attenuation by increasing filter order. 
While nomographs exist for the elliptic and 
Butterworth filters which can be implemented 
by the 2-path structure, a simple approxi- 
mate relationship for the equal ripple 
filter is given in (9), where A(dB) is 
attenuation in dB, Af is transition band- 
width, and N is the total number of all-pass 
segments in the filter. 



A(dB) - (72 Af + 10) N 



2A. HILBERT TRANSFORM FILTERS 



a(n) «=.x(n) + j x(n) 



(10) 



We can cast the 2-path filter as a HT in a 
number of ways but a simple method, akin to 
the transformation performed on a half band 
FIR filter [4 J, is a heterodyne of the half 
band filter to the quarter sample frequency. 
If the original filter response and trans- 
form are denoted by h(n) and H(Z) respec- 
tively, then the heterodyned expressions can 
be easily seen to be those in (11). 

h(n) «> h(n) e'™' 2 - h(n) (j) n (lla) 

H(2) -> H(Z e' ,/? ) -= H(j Z) (lib) 

The transformation of a half band filter to 
an HT filter is accomplished by replacing 
each Z' y of the transfer function with -jZ*- 
Since the polynomials of the all-pass sub- 
filters are second order, this transfor- 



The Hilbert Transform (HT) can be thought of 
as a particular type of half band filter. It 
is modelled and implemented by a wide band 
90* phase shifting network. The HT is often 
used in DSP applications to form the analyt- 
ic signal a(n) as shown in (10), a signal 
with spectra limited to the positive (or 
negative) frequency band. 



mation can be accomplished by changing the 
sign of the coefficient used jn each all- 
pnos sub filtor end by associating a -j with 
the lower path delay. These operations are 
indicated in Fig. 9. 



of Fig. 7. is used for up sampling we would 
perform 3 multiplications per output point. 
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FIGURE 9. HILBERT TRANSFORM ALL- PASS FILTER 
REALIZATION AS TRANSFORMED HALF BAND FILTER 

2B. INTERPOLATION AND DECIKATION FILTERS 

Any half bandwidth filter can be used to 
perform 2-to-l resampling (up or down, 
popularly identified as interpolation and 
decimation) . Finite impulse response (FIR) 
filters can be operated with a polyphase 
partition which does not process inserted 
zero input points (up-sampling) or compute 
the discarded output points (down sampling) . 
Khile this option is not available for the 
general recursive filter, it is an option 
for the recursive all-pass M-path filter. 
The M-path filter can be partitioned into 
polyphase segments because of the interac- 
tion of the internal delays of the M-th 
order subfilters and the delay line of each 
path. This relationship is particularly easy 
to see for the 2-path filter. An impulse 
applied at index n 0 makes a contribution to 
the output at the same time via the upper 
path while an output is not available from 
the lower path till the next index due to 
its extra delay Z'\ During that same next 
index, the upper path offers no contribution 
to the output because there is no Z' 1 path 
associated with its Z' 2 polynomial. Thus the 
filter delivers successive samples of its 
impulse response from alternate paths of the 
filter. 

Down sampling is accomplished by delivering 
alternate input samples to each path of the 
filter operating at the reduced input rate. 
Since the delay of the second path is ac- 
complished by the input commutation the 
physical delay element in that path is re- 
moved. The filtered and down sampled out-put 
is generated by the sum or difference of the 
two , Path responses. Note that the compu- 
tation rate per input sample point is half 
that of the non resampled filter, thus if 
the 2-path filter of Fig. 7 is used for down 
sampling we would perform 3 multiplications 
per input point. 

Up sampling is accomplished by reversing the 
down sampling process. This entails deliver- 
ing the same input to the two paths and corn- 
mutating between their outputs. Both forms 
of resampling are shown in Fig. 10. As in 
the previous example, when the 2-path filter 




FIGURE 10, RESAMPLING 2-PATH ALL- PASS FILTER 

The cascade operation of down sampling and 
up sampling of the complementary outputs 
performs two band maximally decimated quad- 
rature mirror filtering and reconstruction 

Cascading multiple stages of up sampling or 
down sampling filters can achieve high order 
sample rate change with remarkably small 
workload per data point. For instance, a 1- 
to-16 up sampler with 96 dB dynamic range 
for a signal with cutoff frequency of 0.4 f 
requires 2-path filters with successively 
smaller number of all-pass stages of orders 
(3,3), (2,2), (2,1), and (1,1). This re- 
quires a total of 4 2 multiplications for 16 
outputs, for an average workload of approxi- 
mately 2.7 operations per output point. 

2C. ITERATED POLYPHASE ALL-PASS FILTERS 

The all-pass networks described earlier are 
formed by polynomials in Z n and specifically 
for the 2-path filters Z~ 2 . A useful trans- 
formation is to replace each Z" 2 with Z" 2 *, to 
obtain higher order filters exhibiting K- 
fold spectral replication of the original 
filter. For K-2 , the half band filter cen- 
tered at DC becomes a pair of quarter band 
filters centered at DC and f /2. An example 
of spectral replication for the iterated 
filter is shown in Fig. 11. 




FIGURE 11. SPECTRUM OF ITERATED 2-PATH 
FILTER: POLYNOMIALS IN 2* 



Note that this transformation scales the 
prototype filter's spectrum by 1/K, reducing 
both pnnsband width and transition width. 
Thus very steep transition bandwidths and 
narrow filters can be realized by the use of 
low order polynomials of degree 2K obtained 
via delay elements. Energy in the replicated 
spectral regions can be removed by a se- 
quence of filters ! incorporating various 
degrees of resampling and iterated filters 
of reduced degree. 

3. EFFICIENT ALL- PASS FILTER BANKS. 

Efficient spectral partitioning can be ob- 
tained by cascading and resampling the out- 
puts of iterated lovpass and HT filters. For 
instance the spectra of a four channel 
filter bank centered on the four quadrants 
can be formed with a complementary half band 
filter followed by a pair of resampled HT 
filters as shown on the left side of Fig. 
12. The workload for this 70 dB attenuation 
filter set is 3 multiplications per output 
channel . 

A four channel filter bank straddling the 
previous set (ie. centered on the four 
cardinal directions) can be formed with an 
iterated by two complementary half band 
filter followed by both a complementary half 
band resampled filter and a resampled HT 
filter. This spectra of this filter set is 
shown on the right side of Fig". 12. The 
workload for this 70 dB attenuation filter 
set is 2 multiplications per output channel. 
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FIGURE 12. SPECTRA OF FOUR CHANNEL QUADRANT 
CENTERED AND CARDINAL CENTERED FILTER BANK 

4. OTHER ALL- PASS TRANSFORMATIONS 

The prototype complementary all-pass filter 
set can also be transformed to filters of 
arbitrary bandwidth and of arbitrary center 
frequency by standard all-pass transforma- 
tions [6J. The lowpass transformation for 
the half band filter is shown in (12). 



b + Z' 



1 + b Z" 



b « 



l-tan(0 ( /2) 
l+tan(8 0 /2) 



(12) 



This transformation warps the frequency 
variable of the filter so that the -3 dB 



point resides at frequency 0 O . When 6 0 is 
different from r/2 the second order polyno- 
mials of the all-pass subfiltsrs acquire 

coefficients for the Z* 1 ternu The resulting 
filter requires two multiplications per pole 
zero pair which is still half the workload 
of the traditional canonic forms. In addi- 
tion, the pole at the origin (due to the 
delay of the second path) transforms to a 
first order all-pass filter thus converting 
an inactive pole to an active one requiring 
an additional first order stage. Fig. 13 
presents the frequency warped version of the 
filter presented in Fig. 7a. Both this 
filter and its complement are available from 
the warped structure. 




FIGURE 13. LOWPASS TRANSFORMATION OF ALL- 
PASS FILTER 

The bandpass transformation is presented in 
(13). 



-r 1 



1 - c Z 



- c - cos(8 t ) (13) 



This transformation warps the frequency 
variable of the filter so that the center 
frequency resides at frequency 0,. If 0. is 
equal to 0.0 this transformation is equiva- 
lent to the iterated transformation de- 
scribed in 2.C. For any other frequency, 
this transformation converts the second 
order all-pass filters to fourth order all- 
pass structures. This requires four multi- 
plications to form four pole-zero pairs of 
the resultant filter which is still half the 
workload of the traditional canonic forms. 
In addition, the transformation converts the 
real pole (which is an image of the pole 
originally at the origin) to a pair of com- 
plex poles requiring a second order all-pass 
stage. Fig. 14 presents the frequency warped 
version of the filter presented in Fig. 14. 
As earlier, both this filter and its comple- 
ment are available from the warped struc- 
ture. 

Other transformations and geometries [7] can 
be used with the digital all-pass structure 
and the reader is directed to the bibliogra- 
phy of the listed papers for a profusion of 
options . 
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FIGURE 4A. ALL- PASS FILTER STRUCTURE 




(3) 



FIGURE 4B. POLYPHASE ALL- PASS STRUCTURE 
M-l 

2^H n (z") r B . (4) 



H n (Z* K ) 



Y(Z) - X(Z) 



The roots of (3) , the M-th roots of -a and 
its reciprocal, are equally spaced about the 
origin as shown in Fig. 5 for the indicated 
degrees M. Note that ve realize M poles and 
M zeros per multiplication in the all-pass 
stage. This means for instance that a 2-path 
filter offers two poles and two zeros per 
multiplication as opposed to one pole (or 
zero) per multiplication for the standard 
(factored or unfactored) canonic forms. 









FIGURE 5. POLE-ZERO DISTRIBUTION FOR M-TH 
ORDER ALL-PA8S SUBFILTERS 



As ve go around the unit circle we visualize 
that the all-pass subfiltere exhibit rapid 
change in phaee in the neighborhood of each 
pole-zero pair* By proper choice of the pole 
positions, the phase shift of each leg of 
the M path filter can be made to match or to 
differ by multiples of 2w/H over selected 
spectral intervals. This is suggested in 
Fig. 6 for M«2 and H«=4 . The coefficients of 
the all-pass subfilters can be determined by 
standard algorithms [3] or by a new algo- 
rithm developed by harris and d'Oreye and 
reported in a paper recently submitted for 
publication (Signal Processing) . 




FIGURE 6. TYPICAL PHASE RESPONSES OF ALL- 
PASS FILTER PATHS , AND CORRESPONDING KAGNI- 
TUDE RESPONSES. 




FIGURE 7. POLE-ZERO PLOT AND HAGNITUDB 
RESPONSE FOR 2— PATH AND 5 -PATH 
POLYPHASE ALL-PASS NETWORK 
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FIGURE 14. BANDPASS TRANSFORMATION OF 
ALL-PASS FILTER 

5 . CONCLUSIONS 

We have reviewed the form and utility of the 
all-pass polyphase filter structure. We then 
presented a number of all-pass transforma- 
tions which could be easily applied to K- 
path all-pass filter sets. The emphasis has 
been on 2-path networks but the material is 
easily extended to arbitrary M. We alluded 
to a new set of algorithms for designing H- 
path filters which will be reported in an 
upcoming Signal Processing paper. We have 
demonstrated capabilities of this algorithm 
by using it to generate the examples pre- 
sented in this paper. 

The primary advantage of these structures is 
very low workload for a given filtering 
task. Other papers [8,9,10] have discussed 
the low sensitivity to finite arithmetic 
exhibited by these filters. The primary 
impediment to the use of these filters is 
their newness and the lack of available 
design methods for computing filter weights. 
This paper (along with excellent surveys 
listen in the bibliography) address the 
first problem and a forthcoming paper ad- 
dresses the second. 

6. ACKNOWLEDGEMENT 

This work was sponsored in part by the 
Industry/University Cooperative Research 
Center (I/UCRC) on Integrated Circuits and 
Systems (ICAS) at San Diego State University 
(SDSU) and the University of California, San 
Diego (UCSD) . 



7. BIBLIOGRAPHY 

1. Private correspondence between f.j. 
harris and Stan Acks of Hughes Aircraft, 
Radar Systems Group, El Segundo, CA. 1963. 

2. M.R. Schroeder, "Number Theory in Science 
and Communi cations", Springer-Verlag, 1984, 
Chapter 28, Waveforms and Radiation Pat- 
terns, pp. 278-288, 

3. R. A. Valenzuela and A.G. Constantinides, 
"Digital Signal Processing Schemes for 
Efficient Interpolation and Decimation", IEE 
Proceedings, Vol. 130, Pt. G. No. 6, Dec, 
1983, pp. 225-235. 

4. D. Elliot, "Handbook of Digital Signal 
Processing, Engineering Applications", 
Academic Press, 1987, Chapter 3, pp. 227- 
233. 

5. Phillip A. Regalia, Sanjit K. Mitra, and 
P.P. Vaidyanathan, "The Digital All-Pass 
Filter: A Versatile Signal Processing Build- 
ing Block", Proc. of IEEE, Vol. 76, No. 1, 
Jan. 1988, pp. 19-37. 

6. A.G. Constantinides, "Spectral Transfor- 
mations for Digital Filters", Proc. IEE, 
Vol. 117, No. 8, Aug. 1970, pp. 1585 r 1590. 

7. S.K. Mitra, K. Hirano, S. Nishimura, and 
K. Sugahara, "Design of Digital Bandpass- 
/Bandstop Filters with Independent Tuning 
Characteristics", FREQUENZ, 44(1990)3-4, pp. 
117-121. 

8. R. Ansari and B. Liu, "A Class of Low- 
Noise Computationally Efficient Recursive 
Digital Filters with Applications to Sam- 
pling Rate Alterations", IEEE Trans. 
Acoust., Speech, Signal Processing, Vol. 
ASSP-33, No. 1. Feb. 1985, pp. 90-97. 

9. H. Samueli, "A Low-Complexity Multiplier- 
less Half-Band Recursive Digital Filter 
Design", IEEE Trans. Acoust., Speech, Signal 
Processing, Vol. ASSP-37, No. 3. Mar. 1989, 
pp. 442-444. 

10. P.P. Vaidyanathan and Zinnur Doganata, 
"The Role of Lossless Systems in Modern 
Digital Signal Processing: A Tutorial", IEEE 
Trans, on Education, Vol. 32, No. 3, Aug. 
1989, pp. 181-197. 



/************************************** *4************************************^ 
************************** APPe^b*^ ************************** **** 

************************* Least Square Lattice ***************************** 
************************ Noise Cancelling ********* ****************** *^ 

/* Example for ratioroetric approach to noise cancelling */ 
/define LAMBDA 0.95 



void OxiLSL_NC( int 
int 
int 
int 
int 
int 
int 



reset, 

passes, 

*signal_l, 

*signal_2, 

*signal_3, 

*target~l f 

*target_2) 



{ 



int 

static int 
static float 



i f ii f k, m, n, contraction; 
*s_a, *s_b, *s_c 
Delta_sqr, scale, noisejref; 



*out_a, *out_c; 



if ( reset = TRUE) { 
s a « signal_l; 
s_b « signal_2; 
s_c =s signal_3; 
out_a s» target_l; 
out_c « target_2; 
factor « 1.5; 
scale » l.o /4160.0; 

/* noise canceller initialization at time t=0 */ 

nc(0].berr = 0.0; 
nc[0] .Gamma = 1.0; 

for(m«Q; m<NC_CELLS; m++) { 
nc[m].err_a - 0.0; 



Jg. ncfm].err_b « 0.0; 

|^ nc(m).Roh_a = 0.0; 

III; nc[m].Roh_c - 0*0; 

nc(m). Delta « 0.0; 

nc[m] .Fswsqr » 0.00001; 

nc[m].Bsvsqr « 0.00001; 



/ ^mmmmm^mmmmmmmMMmmm^mmmm^mmm END INITIALIZATION 
for(k=0; k<passes; k++) { 
contraction « FALSE; 

for(m=0; m< NC_CELLS; m++) { /* Update delay elements 

nc(m).berrl « nc(m].berr; 
nc[m) .Bswsqrl » nc[m J . Bsvsqr; 

} 



noisejref = factor * log(1.0 - (*s_a) * scale) 

- log(1.0 - (*s_b) * scale) ; 

nc[0].err_a = log (1.0 - (*s_b) * scale); 

nc(0].err_b = log(l.0 - (*s_c) * scale); 



X 



in 



CI 

if! ? 



++s_a ; 
++sjb; 
++s_c; 

ncfO] . ferr 
nc[0] .berr 
nc[0] .Fswsqr 
ncfO] .Bswsqr 



noise_ref ; 
noise_ref ; 

LAMBDA * nc[0]. Fswsqr + noise_ref * noise_ref; 
nc(0] .Fswsqr; 



/* Order Update */ 

for(n=l;( n < NC_CELLS) (contraction 

/* Adaptive Lattice Section */ 

m « n-lj 
ii- n-1; 



FALSE) ; n++) { 



nc[m). Delta *■ 
ncfnj .Delta +«= 
Delta_sqr 

nc[n).fref 
ncfnj.bref = 

ncfnj.ferr = 
ncfnj. berr 



LAMBDA; 

nc[m). berrl * nc[m].ferr / nc(xn] .Gamma 

nc(m] .Delta * ncfmj. Delta; 

-nc(m]. Delta / nc(m] . Bswsqr 1; 

-nc(m). Delta / ncfmj .Fswsqr; 

nc[m].ferr + nc[n].fref * ncf m] .berrl ; 

nc[mj. berrl + ncfnj. bref * ncfmj.ferr; 



ncfn]. Fswsqr «= nc[m]. Fswsqr - Delta_sqr / ncfm] .Bswsqr 1; 
ncfnj. Bswsqr = ncfm] . Bswsqr 1 - Delta_sqr / nc [m J .Fswsqr ; 

iff (nc[n] .Fswsqr + nc(n] .Bswsqr) > 0.00001 | | (n < 5) ) { 

ncfn]. Gamma = ncfm]. Gamma - nc(ra). berrl * nc(m}. berrl / ncfm] .Bswsqrl 
if (nc[n] .Gamma < 0.05) nc[n]. Gamma = 0.05; 
if (ncfnj .Gamma > 1.00) ncfnj. Gamma =_ 1.00; 

/* Joint Process Estimation Section */ 

nc(m].Roh_a *- LAMBDA; 

ncfm].Roh_a +« nc(m].berr * ncfm].err_a / nc[m). Gamma ; 
nc[m).k_a «= nc(m].Roh_a / ncfm] .Bswsqr; 
ncfnj.err_a = nc[m].err_a ~ nc(m].k_a * nc[m).berr; 

nc[m].Roh_c *« LAMBDA; 

ncfmj.Roh_c +«= ncfmj.berr * ncfra].err_b / nc[mj. Gamma ; 

nc[m].k_c «= nc[mj.Roh_c / ncfm) .Bswsqr; 

ncfnj. err_b = nc(mj.err_b - nc[m].k_c * nc[m].berr; 

} 

else { 

contraction = TRUE; 
for(i=n; i<NC CELLS; i++) { 



nc(i] .err_a 
ncfi] .Roh_a 
nc(i j . err_b 
ncfi j .Roh_c 
ncfi j .Delta 
ncf i) .Fswsqr 
ncf i) .Bswsqr 
ncf i} .Bswsqrl 



O.0; 
0.0; 
0.0; 
0.0; 
0.0; 

0. 00001; 
O. 00001; 
0.00001; 



} 



) 



*out_a++ 
*out C++ 



(int) ( (-exp(nc{ii].err_a) +1,0) / scale) 
(int)( (~exp(nc[ii] .err_b) +1.0) / scale) 




Least Square lattice 



|4t 



v. 

Jfck 



APPENDIX C 



/* Normalized Adaptive Noise Canceler */ 



/include <stdlib.h> 
/include <math,h> 

/define TRUE 
/define FALSE 



1 
0 



/define FLO AT 3 2 
/define INT32 
/define VOID 

/define FABSF 
/define SQRTF 



float 
long int 
void 

fabs 
sqrt 



/define MAX(a,b) (a) > (b) ? (a) : (b) 
/define MIN(a,b) (a) < (b) ? (a) : (b) 



/define MIN_VAL 
/define MAX__DEL 
/define MIN_DEL 
/define MAX_RHO 
/define MINJRHO 
fine MIN BSERR 



0.01 

0.999999 
-0.999999 

2.0 
-2.0 

1E-15 



ef struct { 



(ft FLO AT 3 2 
rV| FLO AT 3 2 
•0 FLOAT32 
fR FLO AT 3 2 
Jii FLO AT 3 2 
FLOAT3 2 
^ FLOAT 3 2 
ff FLOAT3 2 
f ! FLO AT 3 2 
W FLOAT3 2 

y float 3 2 

QLANC CELLS; 

HI — 



berr ; 
berr_l ; 
delta; 
err; 
f err ; 
gamma; 
gamma_l ; 
rho; 

delta_l; 
Bserr ; 
Bserr 1; 



typedef struct { 

INT3 2 cc; 
FLOAT 3 2 lambda; 
FLO AT 3 2 min_error ; 

LANC_CELLS * ce 1 1 s ; 
LANC_Context; 



} 



/* number of cells 

/* put in value for lambda 

/* parameter 

/* point to array of ANC CELLS 



extern LANC_ 
LANC_Init( 
INT3 2 
FLO AT 3 2 
FLO AT 3 2 



Context * 



num_cells, 
lambda, 
min__error) ; 



extern VOID 
LANC_Done ( 

LANC_Context 

extern VOID 
LANC_Reset ( 



*c); 



LANC_Context *anc) ; 



/* number of cells 
/* lambda param 
/* min error 



extern FLOAT32 
LANC_Calc( 

LANC_Context *anc, /* input, context handle */ 

FLOAT32 nps, /* input, noise plus signal */ 

FLOAT32 noise); /* input, noise reference */ 

/* The following macros provide efficient access to the lattice */ 

11 



/define 


ANC_CELL 


SIZE 


#def ine 


XBERR 




0 


/define 


XBERR 1 




1 


/define 


XDELTA 




2 


/define 


XDELTA 


1 


3 


/define 


XGAMMA 




4 


/def ine 


XGAMMA 


1 


5 


/define 


XBSERR 




6 


/define 


XBSERR 


1 


7 


/define 


XERR 




8 


/define 


XFERR 




9 


/define 


XRho 


10 



#define berr (* 

#d#€ine P_berr_l (* 

#d<gfine P_berr (* 

#d€|ine berr_l (* 

#dfj|ine Bserr (* 

#de§ine Bserr_l (* 

#df|ine P_Bserr_l (* 

/define P_delta (* 

/define delta (* 

/define delta_l (* 

/ditine P_delta_l (* 

■ «% 

/define err (* 

/define N_err (* 

/define P_ferr (* 

/define ferr (* 

/define gamma (* 

/define P_gamma (* 

/define N_gamma (* 

/define P_gamma_l (* 

/define gamma 1 (* 



p + xBEPvR) ) 

;p + XBERR_1 - ANC_CELL_SIZE) ) 
p + XBERR - ANC_CELL_SIZE) ) 
p + XBERR_1 ) ) 

p + XBSERR) ) 
p + XBSERR_1)) 

p + XBSERR_1 - ANC_CELL_SIZE) ) 

p + XDELTA - ANC_CELL_SIZE) ) 
p 4- XDELTA) ) 
p + XDELTA_1) ) 

p + XDELTA_1 - ANC_CELL_SIZE) ) 
p + XERR)) 

p + XERR + ANC_CELL_SIZE) ) 

p + XFERR - ANC_CELL_SIZE) ) 
p + XFERR) ) 

p + XGAMMA)) 

p + XGAMMA - ANC_CELL_SIZE) ) 

p + XGAMMA + ANC_CELL_SIZE) ) 

p + XGAMMA_1 - ANC_CELL_SIZE) ) 

p + XGAMMA_1) ) 



/define rho (*(p + xRho) ) 



/* 

Name: LANC Init 



Abstract: Create an ANC context 



extern LANC_Context * 
LANG_lnit( 

INT32 num_cells, /* number of cells */ 

FLOAT 3 2 lambda, /* lambda param */ 

FLOAT 3 2 min_error) { /* min error */ 

LANC_Context *anc; /* context */ 

anc = (LANC_Context *)malloc(sizeof (LANC_Context) ) ; 
assert (anc != NULL) ; 
anc->cc = num_cells; 
anc->lambda = lambda; 
anc->min_error = min_error; 

anc->cells = (LANC_CELLS *)malloc(sizeof (LANC_CELLS) * (num cells + 2))- 
assert (anc->cells != NULL); ~ 



} 

/* 



return (anc) ; 



Name : LANCJReset 

Abstract: Reset an ANC context 



exf^rn VOID 
LAfU _Reset ( 

£||LANC_Context *anc) { 

Hi 

>FLOAT3 2 *p; 
J TNT 3 2 m; 

f*t 

JJp = (FLOAT32 *)anc->cells; 
L s for (m = 0; m <= anc->cc; m++) { 
rho = o.O 
err = 0.0; 

ferr = 0.0 
* berr = o.o^ 
berr_l = 0.o| 
delta « o.o] 
delta_l = 0.0^ 
Bserr = anc->min_error; 
BserrjL = anc->min_error ; 
gamma = MIN_VAL; 
gamma_l = MIN_VAL; 
P +- ANC_CELL_SIZE; 

p = (FLOAT32 *)anc->cells; /* Cell # 0 special case */ 

gamma =1.0; " 
gamma_l = i . o ; 



/* 

Name: LANC Done 



Abstract: Delete an ANC context 



*/ 



extern VOID 
LANC_Done ( 

LANC_Context *anc) { 

free(anc->cells) ; 
free(anc) ; 

} 

/* 

Name: LANC_Calc 

Abstract: Calculate 



FLOAT32 
LANC_Calc( 

LANC_Context 

FLOAT32 

FLO AT 3 2 



{ 



*anc, 

nps, 

noise) 



/* input, context handle */ 
/* input, noise plus signal */ 
/* input, noise reference */ 



if 

01 



INT32 

FLOAT32 

FLOAT32 

FLOAT32 

INT32 



m; 

*p; 

B,F,B2,F2; 
qd2 , qd3 ; 
output_cell; 



/* Update time delay elements in cell structure */ 



p = (FLOAT 3 2 *)anc->cells; 

for (m = 0; m <= anc->cc; m++) { 

gamma_l = gamma; 

berr_l = berr; 

Bserr_l = Bserr; 

delta_l = delta; 

p += ANC_CELL_SIZE; 



/* Handle Cell # 0 

P = (FLO AT 3 2 *)anc->cells; 
Bserr = anc->lambda * Bserr_l + noise * noise; 
Bserr = MAX (Bserr, MIN_BSERR) ; 

ferr = noise / SQRTF (Bserr) ; 
ferr = MAX (ferr, MIN_DEL) ; 
ferr = MIN(ferr, MAX_DEL) ; 

berr = ferr; 

rho = anc->lambda * SQRTF (Bserr_l / Bserr) * rho + berr * nps; 
N_err = nps - rho * berr; 

output_cell = anc->cc - l; /* Assume last cell for starter */ 

for (ra = 1; m < anc->cc; m++) { 
p += ANC_CELL SIZE; 



*/ 



B 



= SQRTF (1.0 - P_berr_l * P_berr_l) ; B2 = 1.0/B; 



III 



Hi 



* } 



F = SQRTF(l.o - P_ferr * P_ferr ); F2 = 1.0/F; 

P_delta = P_delta_l * F * B + P berr 1 * P ferr; 
P_delta = MAX(P_delta, MIN_DEL)J 
P_delta = MIN(P_delta, MAX_DEL) ; 
qd3 « l.o - P_delta * P_delta; 
qd2 = 1.0 / SQRTF(qd3); 

ferr = (P_ferr - P_delta * P_berr_l) * qd2 * B2; 
ferr = MAX (ferr, MIN_DEL) ; 
ferr = MIN (ferr, MAX_DEL) ; 

berr = (P_berr_l - P_delta * P_ferr ) * qd2 * F2; 
berr = MAX (berr, MIN_DEL) ; 
berr = MIN(berr, MAX_DEL) ; 

gamma = P_gamma * (1.0 - P_berr * P_berr) ; 
gamma = MAX (gamma, MIN_VAL) ; 
gamma = MIN (gamma, MAX_DEL) ; 

Bserr = P_Bserr_l * qd3; 

Bserr = MAX (Bserr, MIN_BSERR) ; 

rho *= anc->lambda * SQRTF( (Bserr_l / Bserr) * (gamma / gamma 1) ) ; 
rho += berr * err ; — 
rho = MAX (rho, MIN_RH0) ; 
rho = MIN (rho, MAX_RHO) ; 

N_err = err - rho * berr; 
p = (FLOAT 3 2 *)&(anc->cel Is [output cell /* *ANC CELL SIZE 



y return (N_err) ; 



:: : i: 



/* 



QRD 



/include <stdlib.h> 
/include <math.h> 

/define TRUE 
/define FALSE 

/define FLOAT 3 2 
/define INT3 2 
/define VOID 

/define FABSF 
/define SQRTF 



APP"EM)IX D 
*/ 



1 
0 

float 
long int 
void 

fabs 
sqrt 



typedef struct { 

INT32 dummy; 
} LQRD JPEF_CONTEXT ; 

typedef LQRDJPEF_CONTEXT * LQRDJPEF_Handle; 
extern LQRDJPEF Handle 



rt INT32 


NumCells, 


35 FLOAT 3 2 


Lambda , 


&l FLO AT 3 2 


SumErrlnit, 


y| FLO AT 3 2 


Gamslnit, 


Yi FLO AT 3 2 

m 


MinSumErr) ; 


e^tfern VOID 





LQRDJPEF_Done ( 

" LQRDJPEF_Handle hJPE) ; 



exlfefern VOID 
LoMt)JPEF_Reset ( 

II LQRDJPEF_Handle hJPE) ; 



extern FLOAT3 2 

LQRDJPEF_Calc ( 

LQRDJPEF_Handle hJPE, 
FLOAT32 nps, 
FLOAT 3 2 noise) ; 



/* handle 

/* noise plus signal 
/* noise reference 



/define MAX(a,b) (a) > (b) ? (a) : (b) 
/define MIN(a,b) (a) < (b) ? (a) : (b) 



/define LQRD JPEF_CELL_SI ZE 

typedef struct { 

FLO AT 3 2 sinf , sinfJL 
FLOAT32 sinb , sinb_l 

FLO AT 3 2 cosf , cosf_l 
FLOAT32 cosb , cosb 1 



27 



FLOAT32 epsf; 



FLOAT32 epsb , epsb_l ; 



FLOAT32 pief , 
FLOAT32 pieb , 



pief_l 
pieb 1 



FLOAT32 Fserr, Fserr_l, SQRTF_Fserr, SQRTF_Fserr 1; 
FLOAT 3 2 Bserr_l, Bserr_2, SQRTF_Bserr_l , SQRTF_Bserr_2 ; 

FLOAT32 p_l, p_2; 
FLOAT32 gams_l; 
FLOAT32 epsi_l; 



} LQRD JPEF_CELL ; 

typedef struct { 

INT 3 2 

FLOAT32 

FLOAT32 

FLOAT 3 2 

FLOAT 3 2 

FLOAT32 

FLOAT32 

LQRDJPEF_CELL 
} t* LQRD JPEF_Context ; 



NumCells; 
Lambda ; 
SumErrlnit ; 
Gamslnit; 
MinSumErr ; 
SQRTF_Lambda ; /* 
SQRTF_SumErrInit ; / * 
*cells; /* 



/* 
/* 
/* 
/* 
/* 



number of cells 
Lambda 

Initial value for Fserr, Bserr 
Initial value for gams 
Minimum for Fserr, Bserr 
square root of Lambda 
square root of SumErrlnit 
point to array of JPE CELLS 



*/ 
*/ 
*/ 
*/ 
*/ 
*/ 
*/ 
*/ 



I- 



Ir4 



The following macros provide efficient access to the lattice 
Define variables offsets within structure 



#4fifine xSINF 

#<^fine XSINF_1 

/define xSINB 

#dfbbf ine xSINB 1 



#d©eine 
/define 
:me 
:ine 

f W 

/define 
/define 
/define 




xCOSF 
xCOSF_l 
xCOSB 
xCOSB_l 

xEPSF 
XEPSB 
XEPSB 1 



#define xPIEF 
/define xPIEF_l 
/define xPIEB 
/define xPIEB_l 

/define xFSERR 
/define xFSERR_l 
/define xSQRTF_FSERR 
/define xSQRTF FSERR 1 



0 
1 

2 
3 

4 
5 
6 
7 

8 
9 

10 

11 
12 
13 
14 

15 
16 
17 
18 



/define xBSERR_l 19 

/define xBSERR_2 20 

/define xSQRTF_BSERR_l 21 

/define xSQRTF BSERR 2 22 



/define xp_l 
/define xp_2 



23 
24 



/define xGAMSQ_l 
/define xEPSI 1 



25 
26 



/define sinf 
/define sinful 
/define P_sinf 
/define P_sinf_l 

/define sinb 
/define sinb_l 
/define P_sinb 
/define P_sinb_l 

/define cosf 
/define cosfjL 
/define P_cosf 
/define P_cosf_l 

/define cosb 
/define cosb_l 
/define P_cosb 
/define P_cosb 1 

/dBtine epsf 
/dSfine P epsf 

/define epsb 
/define epsbjl 
/define P_epsb 
/dgfine P_epsb_l 

/dgfine pief 

/d^ine pief_jL 

/define P_pief 

/divine P pief 1 

/dfetine pieb 
/define pieb_l 
/define P_pieb 
/define P_pieb_l 

/define Fserr 
/define Fserr_l 
/define P_Fserr 
/define P_Fserr 1 
/define SQRTF_Fserr 
/define SQRTF_Fserr 1 
/define SQRTF_P_Fserr 
/define SQRTF_P_Fserr_l 

/define Bserr_l 
/define Bserr_2 
/define P_Bserr_l 
/define P_Bserr_2 
/define SQRTF_Bserr_l 
/define SQRTF_Bserr_2 
/def ine SQRTF_P_Bserr_l 
/define SQRTF P Bserr 2 



;*(ptr + xSINF)) 
;* (ptr + xsinf_i) ) 

;*(ptr + XSINF - LQRDJPEF_CELL_SIZE) 

;*(ptr + xSINF_l - LQRDJPEF_CELL_SIZE) 

;*(ptr + xSINB) ) 
;*(ptr + xSINB_l)) 

;*(ptr + XSINB - LQRDJPEF_CELL_SIZE) 

;*(ptr + XSINB_1 - LQRD JPEF_CELL_S I Z E ) 

*(ptr + XCOSF)) 
*(ptr + XC0SF_1)) 

*(ptr + XCOSF - LQRDJPEF_CELL_SIZE) 
*(ptr + XC0SF_1 - LQRD JPEF_CELL_SI ZE ) 

*(ptr + xCOSB)) 
*(ptr + XC0SB_1)) 

*(ptr + XCOSB - LQRDJPEF_CELL_SIZE) 
*(ptr + XCOSB_l - LQRD JPEF_CELL_SI ZE ) 

* (ptr + xEPSF) ) 

*(ptr + xEPSF - LQRD JPEF_CELL_S I ZE ) 

*(ptr + xEPSB) ) 

*(ptr + XEPSB_1)) 

*(ptr + xEPSB - LQRD JPEF_CELL_SI ZE ) 

*(ptr + xEPSB_l - LQRD JPEF_CELL_SI ZE ) 

* (ptr + XPIEF) ) 

*(ptr + XPIEF_1) ) 

*(ptr + xPIEF - LQRDJPEF_CELL_SIZE) 

*(ptr + XPIEF_1 - LQRDJPEF_CELL_SIZE) 

*(ptr + xPIEB}) 
*(ptr + xPIEB_l)) 

*(ptr + XPIEB - LQRD JPEF_CELL_SI ZE ) 

*(ptr + XPIEB_1 - LQRD JPEF_CELL_SI ZE ) 



*(ptr + XFSERR) ) 
*(ptr + XFSERR_1)) 

*(ptr + xFSERR - LQRDJPEF_CELL_SIZE) ) 
*(ptr + XFSERR_1 - LQRDJPEF_CELL_SIZE) ) 
*(ptr + XSQRTF_FSERR) ) 
*(ptr + XSQRTF_FSERR_1 ) ) 

*(ptr + xSQRTF_FSERR - LQRDJPEF_CELL_SIZE) ) 
*(ptr + xSQRTF_FSERR_l - LQRDJPEF_CELL_SIZE) ) 

*(ptr + xBSERR_l)) 
*(ptr + XBSERR_2) ) 

*(ptr + XBSERR_1 - LQRDJPEF_CELL_SIZE) ) 
*(ptr + XBSERR_2 - LQRDJPEF_CELL_SIZE) ) 
*(ptr + XSQRTF_BSERR_1) ) 
*(ptr + XSQRTF_BSERR_2 ) ) 

*(ptr + XSQRTF_BSERR_1 - LQRDJPEF_CELL_SIZE) ) 
*(ptr + xSQRTF_BSERR_2 - LQRDJPEF_CELL_SIZE) ) 



/define p_l (*(ptr + xp 1)) 

/define p_2 (*(ptr + xp_2)) 

/define P_p_l (*(ptr + xp_l - LQRDJPEF_CELL_SIZE) ) 

/define P_p_2 (*(ptr + xp_2 - LQRDJPEF_CELL_SIZE) ) 

/define gams_l (*(ptr + xGAMSQ_l) ) 

/define P_gams_l (*(ptr + xGAMSQ_l - LQRDJPEF_CELL_SIZE) ) 

/define epsi_i (*(ptr + xEPSI_l) ) 

/define P_epsi_l (*(ptr + xEPSI_l - LQRDJPEF_CELL_SIZE) ) 

static FLOAT32 
RSQRTF ( 

FLO AT 3 2 X) { 

return 1.0F / SQRTF(x) ; 

/* 

Name : LQRD JPEF_I nit 

Abstract : Create a JPE context 



ejHern LQRDJPEF 


_Handle 


LQlDJPEF Init( 




ffi INT32 


NumCells, 


fU| FLOAT32 


Lambda , 


fiOj FLOAT3 2 


SumErrlnit, 


m FLOAT32 


Gamslnit, 


FLOAT32 


MinSumErr) { 



k<4 

m 



LQRDJPEF_Context *jpe; 

jpe = malloc(sizeof (LQRDJPEF_Context) ) ; 
assert (jpe != NULL) ; 

jpe->NumCells = NumCells; 

jpe->Lambda = Lambda; 

jpe->SumErrInit = SumErrlnit; 

jpe->GamsInit = Gamslnit; 

jpe->MinSumErr = MinSumErr;" 

jpe->SQRTF_Lambda = SQRTF( jpe->Lambda) ; 

jpe->SQRTF_SumErrInit = SQRTF( jpe->SumErrInit) ; 

jpe->cells = ma lloc( si zeof (LQRDJPEF CELL) * (NumCells + 2)) • 
assert (jpe->cells !=NULL); ~ 

LQRDJPEF_Reset ( ( LQRD JPEF_Handle ) jpe) ; 
return ( ( LQRD JPEF_Hand 1 e ) jpe) ; 



/* 

Name : LQRDJPEF_Done 

Abstract : Delete a JPE context 

extern VOID 
LQRDJPEF_Done ( 



} 

/* 



LQRDJPEF_Handle hJPE) { 

LQRDJPEF_Context *jpe = (LQRDJPEF_Context *)hJPE; 

free( jpe->cells) ; 
free(jpe) ; 



Name : LQRDJPEF_Reset 

Abstract : Reset a JPE context 



extern VOID 
LQRDJPEF_Reset ( 

LQRDJPEF Handle 



% A 



Has 



hJPE) { 

LQRDJPEF_Context *jpe = ( LQRD JPEF_Cont ext *)hJPE; 



FLOAT32 
INT32 



*ptr; 
m; 



ptr = (FLO AT 3 2 *)jpe->cells; 

for (m = 0; m <= jpe->NumCells ; m++) { 



sinf 
sinb 
cosf 
cosb 
epsf 
epsb 
pief 
pieb 
P_l 
Fserr 
Bserr_l 
gams l 



0.0F; 
0.0F; 
0. OF; 
0.0F; 
0.0F; 
0. OF; 
0. OF; 
0.0F; 
0.0F; 

jpe->SumErrInit; 
jpe->SumErrInit; 
jpe->GamsInit; 
SQRTF_Fserr = jpe->SQRTF_SumErrInit; 
SQRTF_Bserr_l = jpe->SQRTF_SumErrInit ; 



+= LQRD JPEF_CELL_SI ZE ; 



ptr 
} 

ptr = (FLO AT 3 2 *) jpe->cells; 
gams_l = l.OF; 



/* Cell # o special case 



Name 

Abstract 



LQRDJPEF Calc 



extern FLOAT 3 2 
LQRD JP EF_Ca 1 c ( 

LQRDJPEF_Handle hJPE, 

FLO AT 3 2 nps, /* noise plus signal 

FLOAT 3 2 noise) { /* noise reference 

LQRDJPEF_Context *jpe = ( LQRD JPEF_Cont ext *)hJPE; 



INT32 



m; 



FLO AT 3 2 *ptr; 
FLOAT32 tmp; 

Time update section */ 

ptr = (FLOAT 3 2 *) jpe->cells; 

for (m = 0; m <= jpe->NumCells; m++) { 
/* some of following delay elements are not needed */ 

sinful = sinf ; 

sinb_l = sinb; 

cosf_l = cosf ; 

cosb_l = cosb; 

epsb_l = epsb; 

pief_l = pief ; 

pieb_l = pieb; 

Fserr_jL = Fserr; 
Bserr_2 = BserrJL; 
P_2 = p_l; 

SQRTF_Bserr_2 = SQRTF_Bserr 1; 
SQRTF_Fserr_l = SQRTF_Fserr7 

ptr += LQRDJPEF CELL SIZE; 
} 

/* Order update section */ 
/* Handle Cell #0 */ 

P P tr = (FLOAT 3 2 *) ( jpe->cells) ; /* point to cell # 0 */ 
epsf = noise; 
epsb = noise; 
epsi_l = nps; 



^4 



fi| /* rest of cells */ 

for (m = l; m < jpe->NumCells; m++) { 

ptr += LQRD JPEF_CELL_S I Z E ; /* access next cell */ 

/* Prediction section */ 

P_Bserr_l = jpe->Lambda * P_Bserr 2 + P epsb 1 * p epsb 1- 
P_Bserr_l = MAX(P_Bserr_l, jpe->MTnSumErr) ; ~ ~ ~ ' 
SQRTF_P_Bserr_l = SQRTF(P_Bserr_l) ; 

tmp = RSQRTF (P_Bserr_l) ; /* this comes free on DSP */ 
P_cosb_l = jpe->SQRTF_Lambda * SQRTF_P_Bserr_2 * tmp; 
P_sinb_l = P_epsb_l * tmp; 

tmp = jpe->SQRTF_Lambda * Ppiefl; 

epsf = P_cosb_l * p_epsf - P_sinb_l * tmp; 

P_pief = P_cosb_l * tmp + P_sinb_l * P_epsf; 

gams_l = P_cosb_l * P gams 1; 



P__Fserr = jpe->Lambda * P_Fserr_l + P_epsf * P_epsf; 
P_Fserr = MAX(P_Fserr, jpe->MinSumErr) ; 

SQRTF__P_Fserr = SQRTF(P_Fserr) ; 

tmp = RSQRTF (P__Fserr) ; /* this comes free on DSP */ 

P_cosf = jpe->SQRTF_Lambda * SQRTF_P_Fserr_l * tmp; 

P_sinf = P_epsf * tmp; 

tmp = jpe->SQRTF_Lambda * P_pieb_l; 

epsb = P_cosf * P_epsb_l - P_sinf * tmp; 

P_pieb = P_cosf * tmp + P_sinf * P_epsb_l; 

/* Joint Process Estimation section */ 

tmp = jpe->SQRTF_Lambda * P_p_2 ; 

epsi_l = P_cosb_l * P_epsi_l - P_sinb_l * tmp; 
P_p_l = P_cosb_l * tmp + P_sinb_l * P_epsi_l; 

> 

ptr += LQRDJPEF_CELL_SIZE; /* access next cell */ 

/* Do minimum work for JPE of very last cell 

only four equations are required for prediction section 

PJBserrjL = jpe->Lambda * P_Bserr_2 + P_epsb_l * P_epsb_l; 
PJBserr_l - MAX(P_Bserr_l, jpe->MinSumErr) ; 

SQRTF_PJBserr_l = SQRTF (P_Bserr_l) ; 

tmp « RSQRTF(P_Bserr_l) ; /* this comes free on DSP */ 
P_cosb_l = jpe->SQRTF_Lambda * SQRTF_P_Bserr_2 * tmp; 
P_sinb_l ~ P_epsb_l * tmp; 

/* Joint Process Estimation for last cell */ 
tmp = jpe->SQRTF_Lambda * P_P_2 ; 

epsi_l = P_cosb_l * P_epsi_l - P_sinb_l * tmp; 
P__p_l = P_cosb_l * tmp + P_sinb_l * P_epsi_l; 

gams_l = P_cosb_l * P_cjams_l; 
return (gams_l * epsi_l) ; 



